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Computational Modelling of the Steps 
Involved in Photodynamic Therapy 


L. Therese Bergendahl and Martin J. Paterson 
Department of Chemistry, Heriot Watt University, Edinburgh, 
Scotland 


1. Introduction 


Photodynamic therapy (PDT) is a branch of phototherapy that has seen a surge of interest in 
the last few decades, due to its potential in the treatment of various cancers, infections and 
heart disease.(Bonnett, 2000) This chapter aims to give an overview of the various 
photochemical steps involved in PDT as a cancer therapy, and in particular the challenges 
and insight gained from their theoretical description. After a brief review of PDT in general, 
in a biological and chemical context, the photochemical steps involved will be discussed, 
detailing the computational techniques required to model these chemical pathways 
theoretically. We will detail the methodologies that can currently be applied, as well as their 
limitations of use at present, and areas requiring further development. 


1.1 Review of PDT 


Photodynamic therapy is perhaps best introduced by considering the photochemical 
pathways illustrated in the example Jablonski-type diagram shown below, which describes 
the various processes that would take place within the vicinity of a target cell. 


The photochemical processes involved in PDT starts with irradiation (A). The irradiation 
can take place externally, normally from a laser source, as well as internally, through the use 
of optical fibers, and targets tissue containing a photosensitiser molecule (B). The sensitiser is 
functionalised so that it will accumulate in the cells that are targeted, for example in a cancer 
lesion or an area of diseased skin cells. It is also designed so that it will effectively absorb 
light of a specific wavelength. The absorption process (C) takes place either through the 
absorption of one photon per sensitiser molecule, one-photon absorption (OPA), or, if a light 
source with high fluency is used (e.g., from an intense laser beam), two photons per 
molecule (TPA). In the Jablonski diagram the eigenstates are separated according to their 
energy as well as spin multiplicity, where S denotes a singlet state and T a triplet state. 
These multiplicities, and their relative order, are typical for many organic chromophores 
that are used as photosensitisers, and the transitions taking place between them are 
governed by the various section rules for electronic transitions (e.g. non-parity conserving 
transitions for OPA and parity conserving for TPA). After very fast non-radiative relaxation 
processes, collectively termed internal conversion (IC) (D), the absorption results in the 
generation of the sensitiser in its first excited singlet state, Sı. The sensitiser molecule must 
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Scheme 1. Overview of the photochemical pathways taking place in photodynamic therapy, 
described through the energy flow starting from irradiation of the photosensitiser molecule 
and ultimately ending in cell death. 


then be able to undergo efficient intersystem crossing (ISC) (E) that generates the first 
excited triplet state of the molecule, T1. The sensitiser has to be designed so that ISC takes 
place with high quantum yield. Any radiative processes, where the molecule will return 
to its ground state after emitting a photon, are looked upon as loss processes in PDT and 
should be not effective deactivation channels in an ideal photosensitiser molecule. The 
same applies to any intra-molecular reactive channels in the singlet or triplet manifolds, 
which generate new photoproducts. The comparably long-lived triplet state can react in 
two ways when confined in a biological environment: Type I reactions, which involves the 
creation of radicals via electron or hydrogen transfer reactions, or Type II reactions, where 
an energy transfer process (F) takes place from the triplet state of the photosensitiser, and 
generates the cytotoxic singlet oxygen species (G) from ground state oxygen which is 
always abundant in a biological environment. During the energy transfer process the 
photosensitiser is simultaneously brought back to its singlet ground state where it can, in 
principle, take part in further sensitisation cycles. The singlet oxygen created is 
exceptionally reactive, and is therefore very short-lived in a biological medium where it 
has been shown to be toxic to essentially any biomolecule and organelle within the cell. 
Both pathways ultimately lead to oxidative degradation of the biomolecules within the 
cell, but in most cases Type II mechanisms are dominant, and will be the pathway 
considered in this chapter. Only in environments with low oxygen levels are Type I 
mechanisms observed, which means it will mainly have to be considered for the targeting 
of hypoxic parts of fast-growing tumours. (Bonnett, 2000; MacDonald & Dougherty, 2001; 
Missailidis, 2008) 
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The knowledge of the benefits of sunlight on the general health have been traced back over 
5000 years with references from Egypt and India as well as ancient Greek texts mentioning 
Herodetus and Hippocrates recommending “heliotherapy” to their fellow Greeks and 
documenting the benefits of sunlight on bone growth. (Bonnett, 2000) The development of 
PDT as we know it today has been in progress since the early 1900s, an example being the 
1903 Nobel Prize in Physiology or Medicine awarded to the Danish scientist Niels Rydberg 
Finsen.(Nobelprize.org, 2011) Finsens work was the first of its kind where light was used 
actively and directly as a treatment of disease. The first descriptions of the use of a chemical 
photosensitiser compound in light therapy are also from the early 1900s through the 
observations of Raab as well as Jasionek and Von Trappen.(Ethirajan, Chen et al., 2011) 
Raabs work involved the investigation of bacterial cultures and the effect of an added 
porphyrin dye. He observed completely unexpected cell death in the experiments taking 
place during the daytime, but none during his evening experiments, which lead him to 
draw conclusions about the light activation of his dye molecules and their therapeutic 
outcome. (Allison, Downie et al., 2004) One of the most notable examples from this period is 
the 1913 experiment by the German scientist Meyer, who injected himself with 200 
milligrams of haematoporphyrin, a carbonyl functionalised porphyrin derivative (Fig. 1). After 
sunlight exposure he suffered extreme swelling and experienced photosensitivity for several 
months after the experiment. (Dolphin, Sternberg et al., 1998) A more recent discovery in the 
development of PDT is that a heamatoporphyrin derivative has tumour localising properties 
as well as being phototoxic. This discovery prompted a spur of development between 1960s 
and 1980s which ultimately lead to the first clinically approved drug Photofrin®, based on 
the haematoporphyrin derivative, as well as being the starting point for the development of 
PDT as a promising tool in cancer treatment. (Hasan, Celli et al., 2010) 


The current use of PDT in cancer treatment is typically palliative and for very advanced 
stages of bladder, lung and oesophageal cancers. Its only success as a general treatment is in 
the case of age-related macular degeneration, which is a leading cause of blindness in the 
western world, and this still represents the only area where PDT is the first treatment of 
choice. (Schmidt-Erfurth, 2000) The main limitation in the general application of PDT as an 
anticancer therapy is the poor level of tissue penetration of the light used. Skin cancers are 
the only type that can be treated with relative ease, whilst for lung and oesophagus tumours 
the administration of light normally requires the use of an endoscope. A further problem 
arising if the light dose can be administrated is the targeting of the photosensitiser 
compound. Poor targeting leads to general photosensitivity in patients for weeks after the 
treatment, and even high strength light bulbs can be problematic. General side effects from 
the treatment of internal cancer are nausea and fever, caused by swelling and inflammation 
of the nearby tissue.(Davids & Kleemann, 2011) 


Despite these side effects and even though developments in PDT have not yet resulted in its 
general use as an anti-cancer therapy, it is still an exceptionally attractive goal and research 
continues today along many promising avenues. The attractiveness of PDT as a clinical 
therapy comes from several aspects that could in principle elevate it above current modes of 
treatment. It is a non-invasive therapy that has no toxicity in the absence of a light source, 
and compared to surgery it has excellent outcome with respect to scarring and other 
cosmetic considerations. As the technique generally requires light in the ultraviolet to 
infrared region the radiation damage is minimal, compared to other forms of radiotherapy, 
and it has also been suggested that PDT can trigger and help mobilise the immune system 
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Fig. 1. Haematoporphyrin. The PDT drug Photofrin® consists of oligomers 
haematoporphyrin, of varying size. 


against the tumour.(Brackett & Gollnick, 2011) A further major advantage PDT has, in 
comparison to various forms of chemotherapy, is the inability of the body to develop 
resistance to the treatment as there is not necessarily a specific target biomolecule or 
organelle involved in the mechanism of its cytotoxicity. These facts have stimulated 
continued research as we detail below. One of the most important aspects is a deeper 
understanding of the spectroscopic and mechanistic features as shown in Fig. 1. First 
principles computation has only recently matured to the point where it can be routinely 
applied to some of these problems, and still other challenges remain. Nevertheless such 
studies are now regarded as invaluable in rationalising the existing classes of 
photosensitizers, and can now be used in the design aspect of photosensitisers to have a set 
of desired molecular properties (Arnbjerg, Jimenez-Banzo et al., 2007; Arnbjerg, Paterson et 
al., 2007; Bergendahl & Paterson, 2011; Johnsen, Paterson et al., 2008). 


Any molecule that show promising behaviour for the use as a photosensitiser will have to be 
further optimized extensively in order to consider targeting within the human body, 
toxicity, pharmacokinetics, ease of total synthesis of the molecule, and so on. With 
computational chemistry having reached this stage, systematic studies of molecular 
structure for PDT properties can hopefully aid in the design of drug molecules for 
application in this field. 


2. Computational techniques 


In order to fully explain the fate of an electronically excited species theoretically the 
Jablonski diagram above needs to be expanded, and the energy of each state needs to be 
related to the molecular structure of the system under investigation. This is done by the 
expression of the energy as a function of the degrees of freedom in the molecule. For any 
molecular system, apart from simple diatomic molecules, this will result in a 
multidimensional potential energy surface (PES). This concept, where the motions of the 
atomic nuclei are treated separately from the motion of the electrons, is the result of one of 
the fundamental approximations employed in molecular quantum mechanics: the Born- 
Oppenheimer (BO) Approximation. Conceptually it is based upon the idea that the relatively 
large velocity and low mass of the electrons in comparison to the nuclei means their motion 
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can be decoupled and treated separately, as the dynamics of these species will take place on 
vastly differing timescales. The BO approximation can be used with good results in the 
theoretical description of most areas of chemistry, and the topology of a PES can be used 
routinely to analyse the thermodynamics and kinetics of most thermal reactions. We will see 
however that this approximation fails when describing certain photochemical pathways 
where two PESs can be looked upon as crossing each other, as in IC or ISC when discussing 
PDT, a problem that needs to be addressed for an accurate theoretical description. Thus, 
such treatment of so-called vibronically coupled states is one of the most challenging areas, 
but also one of the most important in Fig. 1 above, as this determines the exchange of the 
initial photo-energy between the electronic and nuclear degrees of freedom. 


2.1 Modelling of the electronic ground state 


When attempting to describe theoretically the PES of the ground state of a large molecular 
system from first principles it is convenient to divide the solution to the electronic 
Schrödinger equation, the total electronic energy, Ea, into various part: 


Ea = Er + Ey + E; + Ex + Ec (1) 


where Er is the kinetic energy of the electrons, Ev is the coulomb interaction between the 
electrons and the nuclei and the coulomb term E; describes the energy the electrons would 
have if they moved independently and if each electron also repelled themselves (i.e., the 
repulsion of the “smeared out” charge clouds of each electron). The final two terms correct 
for the false parts of the latter, Ex, the exchange term, corrects for the fact that electrons with 
the same spin will avoid each other strongly as a result of the Pauli principle, whilst Ec, the 
correlation term, corrects for the smaller correlation of the movement of electrons with 
opposite spin.(Gill, 1998) Most computer codes will calculate an approximate solution to 
these parts of the total electronic energy separately, and the way in which they deal with 
each of them is what distinguishes various methods from each other. That, and of course the 
computational cost of each method (in terms of CPU time). 


2.1.1 Wavefunction methods 


The exact solution to the electronic Schrödinger equation is an intimidating task for any 
system with more than one electron due to the complexity of dealing with interacting 
electrons. The first approximation generally applied is therefore to make a crude 
assumption that the total many-electron wavefunction, ®, can be constructed as a product of 
functions, 6, which describes one electron at the time, i.e. molecular orbitals. This is the 
approach used in the Hartree-Fock method, which is the best approximation available when 
using one-electron molecular orbitals to describe a single electronic configuration, i.e., the 
orbital occupancy of the system. To account for electron spin and indistinguishability the 
orbitals are extended to spin-orbitals and described as a Slater determinant, which satisfies the 
Pauli principle for electrons. The total many-electron Hamiltonian is then used with these 
one-electron orbitals, and the resulting energy is optimised variationally. Crucially however, 
the resulting operator after this treatment is expressed in terms of the molecular orbitals 
obtained as its eigenfunction, creating a non-linear problem that needs to be solved 
iteratively through a so called self-consistent field (GCF) calculation. (Cramer, 2004; Jensen, 
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2007; Szabo & Ostlund, 1996) It should be noted that in the HF method the last term in eqn. 
1 above is completely unaccounted for, and thus this MO method is generally described as 
an uncorrelated method. The electronic configuration described by the HF method need not 
be that of the electronic ground state but if a configuration of a different symmetry is 
constructed then the variational nature of the method can be used to determine the 
optimum molecular orbitals for that configuration. Thus, for an excited state that can be 
qualitatively described by a single electron configuration, of a different symmetry to the 
ground state, then this method is capable of describing such a state. In this so called A-SCF 
approach two problems arise, (i) most excited electronic states are poorly represented by a 
single configuration, and (ii) this method is unable to describe excited states of the same 
symmetry as a lower lying state due to variational collapse. 


One key part in carrying out any computation is the choice of basis functions used to describe 
the one electron orbitals, b, through an expansion of the molecular orbitals via known atom- 
centred functions (ie., the well known LCAO expansion). It is a dilemma to choose 
functions that accurately describe the electron distribution in a given region of space, but at 
the same time are not so numerous as to make the computational cost too steep. Strategies 
that are employed to increase the quality of a given basis set are to include more functions of 
each angular momentum type that are used, as the addition of diffuse functions. (Papajak, 
Zheng et al., 2011) Diffuse functions have a larger radial component than the valence 
functions and can therefore more accurately describe wavefunctions at a larger distance 
from the nucleus, which is useful when considering properties of excited states, as well as 
properties such as polarisabilities, or when dealing with negatively charged species. A 
veritable plethora of basis sets exists for excited state computations. 


The HF method exactly accounts for the fact that electrons with the same spin will avoid 
each other strongly, so-called electron exchange, through the anti-symmetry of the Slater 
determinant. However it completely ignores the smaller correction from the movement of 
electrons with the same spin, this is termed electron correlation. This missing electron 
correlation energy is normally only a fraction of the total energy, but it is often crucial for 
the correct description of a system as it varies greatly with molecular geometry. This is 
especially true for excited states, where there might be more than one electronic 
configuration needed to describe the true bonding situation in the system, so called static, or 
non-dynamic correlation. Thus, static correlation describes the cases where a qualitatively 
correct wavefunction requires several electronic configurations (i.e., a multiconfigurational 
wavefunction). This can either be due to spatial quasi-degeneracies (i.e., small HOMO- 
LUMO gaps) or to spin-symmetry (vide infra). The correlation that takes the instantaneous 
Coloumb repulsion between electron pairs into account, rather than the mean field SCF 
orbital description, is termed dynamic correlation. This can either be recovered from single or 
multiconfigurational reference states. 


The problem of accounting for the missing electron correlation energy is treated differently 
in various methods, but the common aspect of all approaches is the understanding that 
electron-electron interactions can be described more accurately if more than one 
configuration is used to express the final many electron wavefunction, ®. (Jensen, 2007) 
Conceptually the most straightforward solution to the correlation problem is to generate 
more electronic configurations by promoting electrons in occupied orbitals to virtual, 
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unoccupied, orbitals, so that the system is described by a linear combination of all excited 
configurations. This is termed the Configuration Interaction (CI) method: 


Single Double N- fold 
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If the optimisation of the CI expansion coefficients (also known as amplitudes), a;, were to be 
carried out over all possible configurations the CI wavefunction would be the exact 
wavefunction, and hence give an exact solution to the time-independent Schrédinger 
equation within that given one-electron basis set. This is an extraordinarily difficult task as 
the number of possible configuration will grows factorially with system size, and the full CI 
wavefunction normally needs to be truncated in order for it to be routinely used. Even 
though a truncated CI wavefunction will accurately account for a large fraction of the 
missing correlation energy, there is one unfortunate drawback. Truncated CI expansions are 
not size extensive, which means that as the system size gets larger a proportionally smaller 
amount of the correlation energy is returned at a given truncation level. A related issue is 
that the CI method is also not size consistent, ie. the energy is not proportional to the 
number of non-interacting units comprising the system. This is a significant flaw when 
calculating for example dissociation energies, as the energy of two dissociated molecular 
fragments will not be equal to the energy of the two fragments calculated separately, as 
some multiple excitations are omitted.(Piela, 2007) 


A method related to CI is the Coupled Cluster (CC) method. It is arguably the most 
sophisticated and accurate method for retrieving some of the dynamic correlation energy in 


use today. The cluster operator T is related to the exact wave function through a reference 
determinant, usually the HF determinant, through the coupled cluster ansatz; 


vsega EE |g (3) 


where the CC wavefunction is generated through an exponential expansion of the cluster 
operator that generates the substituted configurations from ọ in the same way as in CI 
theory. 


It is assumed that the cluster operator is a sum of excitations operators for all possible 
excitations from occupied to unoccupied orbitals; 


— f 4 

T= +T, tati (4) 
where 7; is an operator for single excitations and so on. If all possible excitations are 
considered the result will be the same as for a full CI treatment, and again this full treatment 
is too computationally expensive to be practical. Yet again a truncation of the wavefunction 
must take place, and this is where the difference to CI becomes evident. Unlike CI methods 
the excitation operators will automatically include important so-called disconnected terms, 
so truncating at T, in the cluster operator automatically includes Car up to 7,” , from the 


Taylor expansion of e”, which results in not only a more accurate description of electron 
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correlation but also generates wavefunctions that are both size extensive and size 


consistent.(Piela, 2007) The disconnected terms like T? can be understood to account for two 
simultaneous pair-wise interactions away from the HF-SCF mean-field description. The one 
disadvantage is that the CC method is not variational (the exponential ansatz unfortunately 
does not generate truncated left and right states in an expectation value of the 
energy).(Crawford & Schaefer, 2000) However this is not as important as the advantages 
gained and therefore CC theory occupies the prime position in correlated wavefunction 
theory. A final advantage the CC approach shares with CI theory is that as the excitation 
operator is expanded one can systematically improve the properties of the molecule 
calculated from the resulting wavefunction (with these methods for ground state properties, 
but see below for excited state properties). However, even with truncation, the 
computational cost of CC computations are very high and the size of the molecules that can 
be investigated are limited, and specialist techniques, such as the use of significant computer 
resources as well as the development of non-iterative versions of the methods, need to be 
employed if larger systems are to be investigated.(Kowalski, Krishnamoorthy et al., 2010) In 
order to routinely be able to probe photosensitiser molecules the size of substituted 
porphyrinic macrocyclic systems different approaches needs also to be considered, for 
example through the use of Density Functional Theory (DFT), although for the core 
chromophores of various photosensitisers coupled cluster theory has set important 
benchmarks against which lower cost methods can be calibrated. 


2.1.2 Density based methods 


In DFT an alternative starting point for the calculation of the total electronic energy is 
employed through the use of the total electron density, p. The main idea of DFT is that the 
total energy can be completely determined by the total electron density of the system 
through a functional that connects it with the energy. (Koch & Holthausen, 2001) The 
advantage to using the density, rather than the molecular wavefunction, is that the problem 
will not depend on the number of electrons present in the system, and therefore avoid 
becoming increasingly complicated as the system increases in size (ie. the energy 
ultimately depends on only 3 spatial coordinates rather than 3n for an n-electron 
wavefunction). In order to develop both an accurate and computationally tractable scheme 
W. Kohn and L.J. Sham developed a fictitious model system where the electrons do not 
interact at all so that they are described exactly by a set of molecular orbitals. Crucially, they 
realised that the density for this non-interacting system can be exactly the same as that of the 
same system with interacting electrons, and this was therefore used to develop the so-called 
Kohn-Sham expression for the ground state electronic energy (where each term is obtained 
from an expression involving molecular orbitals): 


E=ES+E,+E,+ES (5) 


As the expression deals with non-interacting electrons, the kinetic energy term is not 
complete since the exact electronic kinetic energy cannot be written as a functional of the 
density, and its correction is introduced in addition to the exchange and correlation terms, 
which are all then grouped together in the Eycterm. All other terms can be evaluated 
exactly, and so the main problem of DFT is reduced to the estimation of the exchange 
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correlation energy. The way in which various DFT functionals retrieve the exchange 
correlation energy varies, and is what sets them apart. A number of strategies exist and the 
improvement of various functionals is generally derived from both theoretical foundations 
by expanding the parameter set linking the energy with the density (via gradients etc.), and 
the fitting of such parameters to experimental data. 


The Kohn-Sham equation is derived from the Shrédinger equation for the non-interacting 
system and its solution yields the so-called Kohn-Sham orbitals, gi (r), which are defined 
such that: 


AEH © 


These orbitals are then expanded in a basis set, and arranged in a single Slater determinant, 
both of which are completely analogous to HF MO theory. When finding the energy of a 
specific molecular structure using DFT the Kohn-Sham equation is solved self consistently 


using an electron density that initially is a guess in order to calculate Exc, and hence the 
energy of the optimised system.(Koch & Holthausen, 2001) 


2.1.3 Multiconfiguration methods 


The problem with correlation that we have dealt with so far takes the instantaneous 
Coloumb repulsion between opposite spin electrons into account, the dynamic correlation 
(since the Kohn-Sham non-interacting system has a single configuration wavefunction). For 
some systems however there is also a further complication in that more than one electron 
configuration might be needed to fully describe the bonding situation in the system. This is 
especially true for most open-shell systems, of which the oxygen molecule is a prime 
example, as well as in the description of almost all excited states of closed-shell systems. 
This effect is the static or non-dynamic correlation discussed above, and is clearly something 
we need to be able to model if we are to accurately describe all photochemical pathways 
involved in PDT. In the cases where static correlation is important the use of one single 
determinant to describe a many electron wavefunction, with the correct spin, is not always 
possible. There has however been some single-reference methods developed recently that 
can describe some multiconfigurational situations, for example the so-called spin-flip 
approach. In this approach a high-spin triplet reference state is chosen, M, = 1. The target 
state, e.g., the M; = 0 state of an open-shell singlet state, can then be described as a spin- 
flipping excitation from the reference state.(Shao, Head-Gordon et al., 2003) This is a 
promising avenue of attack for the application of DFT to open-shell systems. A more 
common approach in order to deal with static correlation however, is to turn to a strategy of 
taking a linear combination of a number of electron configurations (as in CI theory 
described above), and define a multiconfigurational wavefunction. This wavefunction is then 
optimised variationally by an SCF procedure in which all expansion coefficients are 
variables. That means that the expansion coefficients for the basis functions (the LCAO 
parameters), which describe the molecular orbitals, as well as the configurational expansion 
coefficients will be optimised, which adds flexibility to the wavefunction. Thus the optimum 
orbitals for the given Cl-type wavefunctions are obtained, unlike in standard CI theory 
where the reference orbitals are fixed. The most popular variant of MCSCF theory is the 
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complete active space (CAS) SCF approach. The popularity stems mainly from the fact that the 
choice of which electron configurations that are important to include in the MC 
wavefunction can be difficult, and this is automatically taken care of in CASSCF. In the 
CASSCF approach an active space is selected, that includes specific orbitals, important to 
describe the chemistry one is interested in, e.g., for benzene photochemistry an obvious 
choice is the six m-electrons and six n-orbitals.(Li, Mendive-Tapia et al., 2010) In aromatic 
organic molecules in general, with or without heteroatoms, the frontier molecular orbitals, 
FMOs, (the Highest Occupied Molecular Orbitals, HOMOs, and the Lowest Unoccupied 
Molecular Orbitals, LUMOs) are crucial to describe any photochemical mechanisms, and 
needs therefore to be described accurately. Orbitals of less importance, that are expected to 
stay either fully occupied or un-occupied in all determinants, are kept unchanged, and thus 
no configurations arising from changes to these doubly occupied or empty orbitals are 
included. The remaining orbitals are considered to be active, and a full CI is carried out on 
this subset. CASSCF is exceptionally flexible and can be used to compute essentially any 
part of the PES of excited as well as ground states, including surface crossings, which are 
described later. 


Unoccupied Orbitals 
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Fig. 2. Overview of the CASSCF strategy, where a full CI is carried out in the active space, 
whilst the occupied and unnoccupied orbitals are kept constant in all configurations. 


2.2 Modelling of excited states 


The single reference methods mentioned so far are all constructed for the description of the 
ground state potential energy surface of a molecular system, in particular around regions 
describing stable configurations of the molecule (i.e., PES minima). When dealing with 
photochemical events obviously more than one PES needs to be considered. 


There are many approaches to the theoretical modelling of the properties of a molecular 
system, such as excitation energies, one example being state specific methods. In these a trial 
wavefunction is chosen so that it corresponds to the state we are interested in, which could 
be an excited state rather than the ground state. (Nicolaides, 2011; 2011) This is the idea 
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behind the A-SCF method discussed above. More powerful are methods that can construct 
arbitrary excited state wavefunctions of any symmetry. These are most often obtained as 
higher roots in Cl-type eigenvalue problems, and as such CI methods are the method of 
choice in constructing both qualitative and quantitative excited state wavefunctions. For 
example the simplest qualitative correct excited state wavefunctions can often be 
represented by the CIS approach. Also in this manner CASSCF can generate excited 
wavefunctions of the same accuracy as for ground electronic states. It should be noted that 
such excited states invariably always require such a multiconfigurational wavefunction at 
most points on a photochemical reaction path.(Malmqvist and Roos, 1989; Olsen, 2011) 
Dynamic correlation effects can then be added to such multiconfigurational wavefunctions 
by using the underlying multi-reference configurations for further CI (MR-CI), or 
perturbation theory treatments (e.g., CASPT2).(Pulay, 2011) 


2.2.1 One- and two-photon absorption 


One method for the determination of excited state information that has been very successful 
and popular in recent years is the so-called response function method. The reason for its 
popularity is that it can be adapted to suit DFT, MCSCF as well as CC methods. Thus, most 
molecular properties, such as excitation energies and transition moments can be obtained 
from a response function calculated from a ground state quantity (the wavefunction in 
MCSCF and CC theory, and the density in DFT). In response theory any molecular property 
is expressed in terms of a response function, as a property of a molecule can be defined as 
the response of an expectation value of an operator to a perturbation. (Bast, Ekstrom et al., 
2011) When dealing with perturbations in the form of an external field interacting with a 
molecular system, as when dealing with absorption of light, a time-dependent approach 
needs to be considered. The time-dependent perturbation term is expressed by an 
interaction operator, which is added to the Hamiltonian for the problem. 


The linear response function for operator P with respect to perturbation Q can be written as: 


(v,|P\y,)(v,|Olv.)’ O) 


pS 


An example is the linear response of a system subject to a homogenous electric field of 
frequency @» , where the perturbation operator then will refer to the electric dipole 
operator, p .This equation is identical to the dipole-dipole polarizability tensor at frequency 
@p, i.e. the frequency-dependent polarizability. This has poles (i.e. where the denominator 
is equal to zero) at the excitation energies of the system. The residues of the response 
function (i.e. the numerator at the poles) give the transition moment (related to the 
probability of a transition between states). The poles of the response functions are 
obtained as solutions to an eigenvalue equation, with the eigenvectors giving the 
contribution from orbital excitations to each excited state. This method is very efficient 
since there is no need to construct all the excited states wavefunctions explicitly as in the 
sum-over-states expression in the equation above. All the information regarding the 
excited state is contained in the response function, so there is only a need for one well 
defined wavefunction or density, such as a ground state wavefunction from DFT, MCSCF 
or CC, or density from the Kohn-Sham orbitals in DFT, in order to access information on 
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the excited states. Thus the response function approach allows us to access excited state 
information from nominally ground state methods and take advantage of what they bring 
(i.e., systematic accuracy in CC theory, low cost in DFT) 


As previously mentioned, there is also a possibility of a two-photon absorption process for 
the excitation of the photosensitiser, as per Scheme 2. This has numerous benefits for the 
application to PDT, such as access to the so-called tissue transparency window (i.e., 
avoiding absorbance from water, melanin and haemoglobin) and increased spatial 
resolution for the treatment, thus the accurate modelling of this non-linear absorption 
process is of great importance in the development of ideal photosensitiser molecules. 


Scheme 2. Overview of the TPA process, normally described as taking place through a 
virtual state, v. 


The use of a state specific method for the computation of TPA is inaccurate, as a few state 
model would need to be constructed, so the general strategy is the use of response theory 
through higher order terms of the response function. For TPA the quadratic response 
function needs to be considered under two dipole perturbations, one of frequency œs and 
one of frequency @. The energy representation of the generalised quadratic response 
function is then described as: 


Ply. yw, Blw (wv, Rw) 


v,jži (a, z. O, + o)(a, = o,) 


(8) 


When the operators P, Q and R are dipole operators then this expression is identical to the 
first dipole hyperpolarisability tensor at frequencies œp and œc. This frequency-dependent 
hyperpolarisability can describe non-linear phenomena involving two photons of 
frequencies ©, and @, respectively, as well as TPA information, where two identical photons 
are absorbed simultaneously. (Paterson, Christiansen et al., 2006) The residues of the 
quadratic response function can be shown to give the conventional expression for the two- 
photon absorption transition matrix element, M, between states i and j through a virtual state 
(written as linear combination of all real eigenstates v): 
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where o the frequency of the irradiating light, and @, is the excitation frequency of state v. 
The labels a and 6 refers to the Cartesian coordinates x, y and z of the TPA tensor. The 
transition matrix element is then used to calculate the rotationally averaged two-photon 
absorption cross section, 6. 


Even though the formal energy expressions for the quadratic response function appears to 
involve a tedious sum over all states of the molecule to describe the virtual state v above, as 
mentioned previously the actual computational machinery of response theory calculates the 
properties directly without constructing the states explicitly which saves both time and 
computational effort. Additionally, one does not need to make any a priori estimate of the 
importance of certain states to the TPA tensor. Not only that, this technique has also been 
proven to be exceptionally accurate in the prediction of TPA spectra. When response theory 
is used with a CC ground state wavefunction quantitative agreement with experimental 
results can be achieved, and for a DFT ground state density very good agreement is also 
often obtained.(Arnbjerg, Jimenez-Banzo et al., 2007; Hattig, Christiansen et al., 1998) When 
it comes to DFT however, the choice of functional is crucial for the application of quadratic 
response theory, as it has to be able to model excitations to and from the virtual state, which 
essentially means it has to be considerably more robust than the functionals routinely used 
together with linear response theory for one-photon absorption. If this requirement is 
satisfied quadratic response DFT has been shown to be very accurate, and even comparable 
to high-level CC calculations for TPA properties.(Arnbjerg, Paterson et al., 2007; Johnsen, 
Paterson et al., 2008; Paterson, M. J., Christiansen et al., 2006) This is very encouraging, 
considering the size of the photosensitisers we are interested in modelling for the 
description of TPA in PDT (vide infra). The improvement of non-linear absorbing qualities of 
photosensitiser molecules can involve many complicated synthetic routes(Drobizhev, 
Karotki et al., 2002; Drobizhev, Meng et al., 2006) As the benefits of the use of TPA in PDT 
are plentiful, a theoretical input will be invaluable, as it can give important information on 
molecular structure-TPA relationships which are crucial for the development of TPA PDT 
on a daily basis. This has for example been demonstrated in an experimental and theoretical 
study, where the porphycene macrocycle was shown to be a promising photosensitiser 
system for TPA PDT.(Arnbjerg, Jimenez-Banzo et al., 2007) 


We have carried out further computational TPA studies of this porphycene molecule and 
shown that the non-linear absorption can be effected to a somewhat unexpected degree by 
the substitution of heteroatoms in the macrocyclic core.(Bergendahl & Paterson, 2011) The 
structural changes afforded very effective absorption in the so-called Soret-region of the 
TPA spectra (an unusual feature for the porphyrin photosensitisers so far investigated for 
TPA PDT) without having a large impact on the one-photon absorption spectra. We 
investigated a set of dioxaporphycenes, where the pyrrolic nitrogen atoms were substituted 
for oxygen atoms systematically as per Scheme 3. 


After optimising the geometries of the molecules, the OPA and TPA was calculated using 
linear and quadratic DFT response theory respectively, with the CAM-B3LYP functional and 


www.intechopen.com 


174 Molecular Photochemistry — Various Aspects 


Pc-O1 Pc-O2 Pc-O3 


%g@ 2° 2%, ¢@ og” 26 oa” 
d a 29” a aa a" “a. 
3 ə ə ə ə 3 
Can C,, Ca 


Scheme 3. The optimised geometries of the three dioxaporphycenes investigated, with 
resulting symmetries. N atoms in blue, and O atoms in red. 


6-31G* Pople-type basis set. Previous work has shown that the CAM-B3LYP functional is 
very accurate at modelling TPA, as it has correct asymptotic behaviour in the exchange- 
correlation functional, which gives a good description of absorption to and from the virtual 
state.(Paterson, Christiansen et al., 2006; Rudberg, Salek et al., 2005) 
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Fig. 3. Summary of the OPA spectra for the dioxaporphycenes investigated. Note the 
absorption band at 300 - 350 nm, the Soret-band, and at 520 - 580 nm, the Q-band, which are 
characteristic bands for the porphyrin family 


The OPA results are summarised in the simulated absorption spectra in Figure 3, and show 
clear Soret- (~ 300-350 nm) and Q-bands (~ 500 - 600 nm) for all species. As expected, there 
are no major effects observed on the spectral profile upon substitution of the core 
heteroatoms, compared with porphycene. The TPA spectra are presented in Figure 4. It is 
important to note that the absorption wavelength indicated is the total for the two-photon 
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process, i.e. the actual excitation takes place at twice the wavelength. The spectra all show a 
main absorption in the Soret-region (~ 250 - 270 nm), and there is a significant disparity 
between the calculated TPA cross-sections for the four isomers. These two facts can be 
explained in terms of one-photon resonance enhancement of the absorption. This 
phenomenon is best discussed in terms of the sum-over-states expression in the equation for 
the transition matrix element M, as seen previously, where the magnitude of M depends on 
the difference between the excitation frequency of the v-th state and the frequency of the 
total excitation, from the magnitude of the denominator in the expression. As the virtual 
state can be described as a linear combination of all real eigenstates in the system, a TPA via 
a virtual state that can be described mainly in terms of a real state will be significantly 
enhanced. Thus, an allowed OPA Q-band transition accidentally degenerate at half the Soret 
excitation energy will greatly enhance the Soret TPA transition. These results suggest that 
the tuning of TPA does not have to be carried out through large structural changes, and 
small electronic factors can be used to improve the magnitude of the TPA cross-section, 
hence broadening the scope for the development of TPA photosensitisers in PDT. 
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Fig. 4. Summary of the TPA spectra. Note the change in scale between the two TPA cross- 
section axes, 


2.2.3 Photochemistry on excited state potential surfaces 


Even though the information we gain from response theory can accurately predict the 
electronic spectra for both one and two photon absorption processes, it is insufficient to give 
a full overview of the whole reaction path. This is true when we are modelling reactive 
photochemistry in general, and the paths involved in PDT in particular where IC and ISC 


www.intechopen.com 


176 Molecular Photochemistry — Various Aspects 


plays a significant role in the understanding of the reaction path. The modelling of reaction 
paths which cross multiple surfaces have been studied both experimentally and 
theoretically for number of years, and the understanding of the pathways involved in the 
photochemistry of small molecules have expanded in great detail.(Zewail, 2000; 2000) 


The non-radiative IC pathways have been known to exist for a very long time, an early 
description being Kashas rule, which states that radiative transitions and ISC will take place 
from the lowest energy excited level of that spin multiplicity due to the very short time-scale 
of the IC processes between upper electronic levels.(Kasha, 1950) More recently however, 
the description and understanding of these events have experienced a dramatic increase 
through the development of electronic structure methods to theoretically model these 
events. It is known that IC events are normally dominated by so called conical intersections, 
where two, or more, surfaces cross each other. Conical intersections are mathematical 
objects that have been known since the 1930s, and many excellent reviews exist on the topic. 
(Matsika & Krause, 2011; Worth & Cederbaum, 2004; Yarkony, 1998) They were well known 
early on to exist in high symmetry species, and manifest through phenomena like the Jahn- 
Teller effect, showing that molecules in orbitally degenerate states undergo non-totally 
symmetric distortions to remove the degeneracy.(Paterson, Bearpark et al., 2005) However, 
advances in electronic structure theory in the last 20 years or so have highlighted that 
accidental conical intersections, not requiring high symmetry, are exceptionally common in 
polyatomic systems, and especially in heteroaromatic molecules, which is of interest from a 
photosensitiser point of view. 


Potential Energy A 
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Fig. 5. When the energies of the crossing surfaces are plotted agains two coordinates making 
up the so-called branching space, the surfaces meet in a characteristic double cone shape 
forming a funnel along the reaction path. 
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The experimental evidence for the presence of conical intersections often require advanced 
set-ups and generate complex data sets, such that a theoretical input is often required for 
full characterisation.(Hadden, Wells et al., 2011; Oliver, King et al., 2011) In non-radiative 
processes, where a system moves from potential energy surface to another, this population 
transfer is most efficient at points of degeneracy (i.e., conical intersection seams). 
Furthermore, the timescale for such transitions is frequently found to be ultrafast (i.e., 
femtosecond) as a single molecular vibration can cause a change in electronic state. Clearly 
non-adiabatic terms are very important, as the dynamical processes involve the coupling of 
nuclear and electronic motion to be taken into account if the full reaction path is to be 
understood.(Worth & Cederbaum, 2004) 


Modern computational codes can straightforwardly locate stationary points on a surface, 
but in the case of two surfaces crossing, the actual crossing point is not a stationary point 
of either of the surfaces. The intersection seams are of dimension 3N-8 for an N-atom 
molecule (every point on the seam being a conical intersection geometry). Modern 
optimisation algorithms optimise to points of zero gradient in this so-called intersection 
space, generating minimum energy crossing points.(Paterson, Bearpark et al., 2004; 
Yarkony, 2005) These act as photochemical funnels for the reaction in the upper electronic 
state. It is important to note however that this minimum energy crossing point might not 
be the only one of importance when it comes to the generalisation of various reaction 
paths, but a good characterisation of the full seem is crucial.(Paterson, Robb et al., 2005) 
The procedure that is needed to find the minimum energy crossing point is actually fairly 
straight forward, at least in principle, and involves computing both the energy and its 
gradient on both surfaces, as well as the non-adiabatic coupling, which essentially 
describes the extent to which the electronic states are coupled by nuclear motion (i.e., a 
break-down of the BO approximation discussed earlier). Multiconfigurational methods 
are the most common, and appropriate, for the calculations of the IC pathways, as they 
can handle the multi state character needed near a crossing point. However, as the 
ubiquitous nature of conical intersections, and their importance in the description of 
photochemistry, in molecules with biological importance, such as photosensitiser 
molecules, becomes more evident there is a need to develop methods that can handle 
systems of a larger size. The use of TD-DFT, which is successful for excited states in the 
vertical region, in the characterisation of crossing points have been investigated, and has 
been shown to be of some use for the location of minimum energy crossing points 
between excited levels (i.e., not with the ground electronic state) They can not yet describe 
excited state surfaces accurately in general, so care must be taken in their application to 
problems in reactive photochemistry as opposed to spectroscopy. (Dreuw & Head- 
Gordon, 2005) 


When it comes to the description of surface crossings of larger molecules, alternative 
multiconfigurational approaches are required. This can be carried out in two general ways: 
the truncation of the computational method we want to use, or the truncation of the 
molecule to a representative model system. Bearpark et al, investigated a collection of 
photochemical reactions involving multi-state processes using a CASSCF starting 
point.(Bearpark, Ogliaro et al., 2007) They primarily considered so called hybrid methods 
where the molecular system is partitioned so that one area, such as the reaction centre or 
chromophore, is treated with a more accurate, and computationally more expensive, model. 
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In addition they investigated a conventional method, where the whole system is treated at 
the same level of theory, although biased towards a particular region of the molecule. At the 
moment the practical limit of the use of CASSCF lies in the region of 14 active electrons in 14 
active orbitals. Methods that can combine a quantum mechanical method (QM) with a 
Molecular Mechanics method (MM) (which determines molecular properties through the 
use of empirically built force fields) are generally referred to as hybrid methods. Hybrid 
methods, such as ONIOM investigated in this study, rely on the fact that some parts of the 
molecule are more crucial in the representation of the chemical reaction we want to study, 
whilst the remainder are viewed as playing a minor role. The ONIOM version of QM/MM 
method has the advantage of being very general, as in principle more layers of levels of 
theory can be applied. This method was found to be appropriate for excited state potential 
energy surfaces, provided the chromophore could be accurately separated from the full 
molecule. This highlights the importance of these models in the investigations of 
chromophores in for example a protein or solvent environment. When it comes to molecules 
where the chromophore is essentially the entire system however, these methods were not 
always appropriate. This is a key factor when dealing with photosensitisers for use in PDT, 
as they normally are of a large structure incorporating many aromatic features, such as 
porphyrins, where the whole system is crucial for the modelling of its photochemistry. For 
molecules such as this the strategy would be to turn to models where the method itself, 
rather than the molecular structure, has been truncated, such as the RASSCF 
approach.(Merchan, Orti et al., 1994; 1994) The RASSCF method is analogous to CASSCF 
where the active space has been divided into smaller parts, as to reduce the total number of 
electron configurations generated. The active space is separated according to maximum and 
minimum occupancy criteria, where the RAS 1 space only allows a limited number of 
vacancies, holes, and the RAS 3 space only allows a limited number of electrons. In the 
remaining RAS 2 space a full CI is carried out, as in the general CASSCF approach. In this 
study it was found that the usefulness of the various methods studied was very system- 
dependent. However, providing the partitioning of the system, or selections of orbitals in 
the RASSCF case, were adequate the calculations lead to small errors when modelling the 
excited state surfaces, when compared to CASSCF, on a 14-electron system.(Bearpark, 
Ogliaro et al., 2007) 


We are currently investigating the use of RASSCF for analysis of surface crossings in 
molecules of the size of common photosensitiser molecules. The RAS space chosen here 
highlights the difficulty in the use of multiconfigurational methods for large molecules, as 
essentially the whole molecule forms a crucial part of the chromophore (Fig 7). The selection 
of the orbitals to be included in RAS-2 has to be founded on photochemical background of 
the problem. In the case of porphyrin derivatives it is well known that the two dominating 
absorption bands in the spectrum (the Soret- and Q-bands) arise from excitation from the 
HOMO and the HOMO -1 to the LUMO and the LUMO+1. The use of these FMOs to 
describe the absorption characteristics of porphyrin systems is termed the Gouterman Four 
Orbital (GFO) model and it is obvious that an accurate description (i.e., full CI in the RAS 2 
space) of these orbitals will be crucial for an accurate RASSCF calculation (Fig 7). A 
methodical description of a 26-electron 24-orbital system like these is already a challenging 
task, and in order for these methods to be useful for the investigations into photosensitiser 
molecules with various specialist properties, like large TPA cross section, advances in both 
electronic structure codes and computational performance will be crucial and lead to 
important insight for these types of investigations. 
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Fig. 6. Overview of the RASSCF strategy, where the active space is divided into three parts, 
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Fig. 7. Overview of the orbitals chosen in the description of the porphyrin chromophore, 
where the whole n-system have to be accounted for. RAS 2 includes the frontier molecular 
orbitals, FMOs, of the porphyrin macrocycle, as these have been shown to be vital for an 
accurate description of the absorption in molecules of this type. 
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In the computation of the crossing of states with different spin symmetry, ISC, it is necessary 
to consider spin-orbit coupling, the coupling of the net spin and orbital angular momenta of 
the electrons. This can be viewed as a small amount of spin being mixed between two states, 
which lead to intensity from ‘spin-allowed’ transitions getting transferred to formally ‘spin- 
forbidden’ transitions, hence facilitating the increased probability of the latter. The process 
of ISC can be modelled through multiconfigurational methods, as described above for an IC 
process. It should be noted that for weak spin-orbit coupling in photosensitisers the 
geometries at which the singlet and triplet surfaces cross are the dominant features that 
need to be understood. These crossings are 3N-7 dimensional seams in contrast to conical 
intersections due to the non-adiabatic coupling being zero by symmetry. Standard conical 
intersection algorithms can be routinely applied to find minimum energy points on these 
ISC seams. The use of RASSCF to model this is now possible given recent advances in the 
RASSCF method, in particular the use of gradient driven RASSCF. 


There has also been considerable development in the use of DFT to model a change in spin 
in the ground state of a molecule, mainly spurred from the interest in the modelling of 
chemical reactions involving a change of spin in transition metal compounds, which can be 
of interest from the point of view of photosensitiser molecules.(Poli & Harvey, 2003) Even 
though DFT has had most success as a ground state method, the modelling of crossings 
between spin-states can be carried out effectively, as the event can be viewed as a single 
surface reaction. This is done by considering the reactants and products of the 
photochemical reaction as being located on two diabatic surfaces that have significant 
mixing in the area of the transition state. If the surfaces of the two spin-states are described 
accurately by DFT, the method can be used with good results to predict the energy of the 
crossing point between them. Of course this only applies to the crossings between the lowest 
energy states of each multiplicity. The main issue with DFT however is still the choice of 
functional to apply to the problem, as there is no direct recipe even for what class of 
functional will preform better in these cases. (Reiher, Salomon et al., 2001) 
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Fig. 8. Two diabatic surfaces, each desribing a spin state, crossing in an area which can be 
viewed as the adiabatic TS. 
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3. Photosensitisers in PDT 


As in the development of any drug molecule there are many aspects to consider in the 
advancement of an ultimate photosensitiser molecule, such as pharmacokinetics, targeting 
and cytotoxicity. However, in the case of a photosensitiser molecule, there is also the 
absorption qualities to consider and, in particular, the location of the wavelength of 
maximum absorption in relation to the of tissue transparency window. 
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Fig. 9. The Porphyrin free base, with crucial structural sites indicated. 


All modern photosensitisers used in practice today are based on the porphyrin macrocyclic 
structure, and PDT with porphyrin derivatives have been under development since the 1960s 
(Fig 9).(Peng & Moan, 2003) In the body, the most recognisable porphyrin derivative that exists 
naturally is coordinated to Fe in the haem structure, present in the haemoglobin protein. The 
unchelated porphyrin structure is known as protoporphyrin IX (PpIX), which has been shown 
to be an effective photosensitiser.(Kennedy, Pottier et al., 1990) The strategy for PpIX treatment 
is to administer a synthetic version of 5-Aminolaevulic acid, the biosynthetic precursor to 
PpIX, which leads to a build-up of PpIX in the targeted area. In Europe, PpIX is approved for 
use in PDT of actinic keratosis and basal cell carcinoma.(Morton, McKenna et al., 2008) The 
most common photosensitiser in use today, however, is still the heamatoporphyrin derivative 
(HpD), which consists of a mixture of oligomeric porphyrin chains. In its purified form it is 
known as Photofrin®, which was the first PDT sensitiser to be approved for clinical 
use.(Elsaie, Choudhary et al., 2009) Photofrin or HpD is approved in Europe, the USA, Canada 
and Japan and is used mainly for cervical, bladder and gastric cancers. 


Although Photofrin®, and other porphyrin derivatives, show excellent singlet oxygen 
quantum yields the longest wavelength at which they are excited is at 630 nm.(Dolphin, 
Sternberg et al., 1998) This maximum absorption is within the transparency window of 
human tissue as desired but not at an optimum position. Light of a longer wavelength 
would penetrate deeper and the clinical use of the photosensitiser would be expanded. This 
has lead to extended research in ways to manipulate the wavelength of maximum 
absorption in porphyrin related molecules and has lead to substituted porphyrins being the 
most studied class of sensitiser molecules.(Peng & Moan, 2003) One strategy to alter the 
porphyrin structure is to add substituents in the meso-position. Meso-substituted 
porphyrins were investigated by Oliviera and co-workers, and they could demonstrate that 
unsymmetrical substitutions on the porphyrin molecule had neglible effect on the 
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absorbance as well as on the fluorescence spectrum. The effect on the singlet oxygen 
quantum yield was also neglible. These were significant results, as they points out a 
porphyrinic site that is available for potential modification with drug targeting in mind, 
without necessarily affecting the chromophore itself.(Oliveira, Licsandru et al., 2009) 
However, most known porphyrins do not absorb at longer wavelengths than Photofrin, so it 
is likely that interest as well as importance of pure porphyrins for use a photosensitiser 
agents in PDT will decrease in the future, whilst concentration is being shifted towards 
porphyrin derivatives, so called second generation photosensitisers. 


Substituting one N atom in the porphyrin core for a heavier S or Se atom, in order to facilitate 
larger spin-orbit coupling, have little effect on the absorption spectrum, apart from a shift in 
the absorption maximum to slightly longer wavelengths (~665 nm), (Rath, Sankar et al., 2005) 
(Detty, Hilmey et al., 2002) however without any effect in the generation of singlet oxygen. 
Substitution of a further N atom with S or Se moves the absorption maximum to longer 
wavelengths, but shows a negative effect in the generation of singlet oxygen.(Detty, Gibson et 
al., 2004) Even though the singlet oxygen generation is decreased in the structures where two 
core heteroatoms have been substituted, they have been proven to be more effective than 
meso-substituted pure porphyrin cores in in vivo studies.(Detty, Gibson et al., 2004) The 
avenue for the use of core-substituted macrocycles in TPA PDT is still in its early stages of 
investigation, but does show some promise.(Bergendahl & Paterson, 2011) 


Porphyrin isomers that has one of their B-carbons and pyrrolic nitrogen atoms switched are 
termed confused porphyrins and they show a maximum absorption at longer wavelengths 
than porphyrin (~ 730 nm).(Furuta, Asano et al., 1994) Even though these types of compounds 
have potential as photosensitisers, with moderate singlet oxygen yield, their difficult total 
synthesis is a major obstacle for their development as photosensitiser agents in PDT. 


A molecular structure related to the porphyrin core can be found in the porphycene 
compounds, a further reduced porphyrin isomer. The porphycene core structure have been 
shown to experience a large absorption at slightly longer wavelengths than porphyrins (~ 
633 nm) which is thought to be due to lower molecular symmetry, and in vitro studies prove 
a modest singlet oxygen quantum yield.(Stockert, Canete et al., 2007) Porphycenes have also 
been shown in in vivo studies to exhibit selective uptake in membrane related organelles, 
which was proven to ultimately lead to apoptosis of the cells after application of radiation. 
(Stockert, Canete et al., 2007). The importance of the porphycene system has been suggested 
to lie mainly in the application of TPA PDT, as early theoretical studies suggests its 
superiority at effective TPA compared to many related porphyrin derivatives. 


As the main focus of the investigations into effective photosensitiser system had been 
dominated by synthetic developments of porphyrin related systems, this is also the avenue 
that recent theoretical efforts have gone down. Hence the theoretical methods mentioned in 
this chapter have been discussed with systems like these in mind. It is however important to 
keep in mind that non-porphyrin photosensitisers also are being investigated, including 
various dye molecules and systems on the nano-scale such as semiconductor quantum dots 
and fullerenes.(Bakalova, Ohba et al., 2004; Mroz, Tegos et al., 2007) As the field of 
theoretical chemistry develops, with improvements in algorithms as well as advances in 
computational hardware, the input into PDT research has potential of becoming very 
important in the future. 
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4. Oxygen in PDT 


The species we have so far referred to as ‘singlet oxygen’ in our discussion on PDT is the 
first excited state of molecular oxygen (a!Ag). Reactions involving interactions from this state 
play a key role in many photodegradation processes in photochemistry as well as 
photobiology and various applications exists, aside from PDT, examples including 
mechanistic cell death studies, bleaching processes and general organic synthesis 
applications. (Schweitzer & Schmidt, 2003) 
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Fig. 10. Molecular orbital diagram for oxygen. The occupation of the antibonding n-orbitals, 
which give rise to the various excited states of Oz is indicated. 
Adapted from (Paterson, M. J., Christiansen et al., 2006) 
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The a Ay -X `E, transition is arguably the most forbidden electronic transition in nature, 
with violation of not only the spin selection rule, but parity and orbital selection rules as 
well, and in fact all the radiative processes of the singlet states of the O2 molecule are 
forbidden to a various degree. This leads to an extraordinary metastability of these states in 


the gas phase, with lifetimes of 45 min for the 1A, state and 7.1 sec for the 'Z; state.(Kasha & 
Brabham, 1979) Crucially, however, the degree that these radiative transitions are forbidden 
decreases rapidly as soon the molecule is perturbed by other atoms or molecules. The nature 
of this perturbation is generally thought of as a weakly bound collision complex, M-O2, 
where M could be a number of species, including solvents or O2 itself. In the case of 
sensitisation O2(X°2,) forms a collision complex with an excited state of a photosensitiser 
molecule, as will be described more later. The weakness of this interaction presents 
problems for computational modelling. Firstly, many possible conformations may be 
possible and some dynamical sampling will be required. There are at present no good force- 
fields to model this with MD. Secondly, methods such as DFT do not generally describe 
such non-covalent interactions very well, although empirically corrected functionals 
describing such dispersion effects have recently been developed. 


As illustrated in the molecular orbital diagram in Figure 10, 14 electrons (including the 1s? 
pairs) of the O2 molecule can be paired up in the lowest lying molecular orbitals. This leaves 
two remaining electrons that have to occupy a doubly degenerate set of orbitals (xg). There 
are several possible ways to accomplish this. This is a basic and simple view of the electronic 
structure of Oz, but it highlights some important aspects of the system, such as the fact that 
the differences in the three lowest lying states in O2 are only a result of the pattern in which 
the final two electrons occupy two orthogonal and degenerate molecular orbitals. It also 


highlights the fact that the photochemically important 4'A,and X°2, states are both open- 
shell states, i..e, states involving partially occupied orbitals. This open-shell character is 
crucial to consider if we are to attempt to accurately model Oz, and M-Oy collision 
complexes computationally.(Schmidt & Schweitzer, 2003) (Paterson, Christiansen et al., 
2006) When modelling the oxygen molecule, and its role in PDT, the important states are in 
fact all of open-shell character to some extent, and often require several electronic 
configurations to describe them. This open-shell character gives rise to physical phenomena 
that a computational method needs to be able to account for. One example is so called spin- 
polarisation, which accounts for the fact that up and down spin electrons will respond 
differently to an excess of spin density, which in turn affects the spatial distribution of up 
and down spin electrons, meaning that it is no longer equal for each as is the case for a 
closed shell wavefunction. The wavefunction should also be constructed so that spin 
angular momentum is accounted for properly, i.e., the wavefunction should be an 
eigenfunction of the spin angular momentum operators describing the total spin (S) and the 


projection of spin along one direction (the operators $ and S, with eigenvalues S(S+1) and 
Ms respectively). This ensures that the wavefunction describes a pure spin state, without 
being contaminated by terms from other spin states that describe additional unpaired 
electrons in the model, a problem referred to as spin-contamination. (Jensen, 2007) High- 
spin open-shell states such as *X `E, can be treated qualitatively using a single- 
configurational model, while low-spin open-shell states (including a'A,) cannot. If one 
wishes to use a single configurational model then both the polarisation and contamination 
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cannot be simultaneously accounted for. One can either use a restricted approach in which 
spin up and spin down electrons are maximally paired and share a common spatial 
component, or one can use an unrestricted approach in which each up and down spin 
electron have different spatial orbitals. In the restricted approach the single configuration is 


automatically an eigenfunction of S, with an eigenvalue corresponding to half the excess 
number of spin up electrons (using the convention of the open-shell electrons having spin 
up), and the eigenvalue of §?, corresponding to a pure value of S, given by the sum of spin 
up quantum numbers for all open-shell electrons. Since these open-shell electrons see the 
same spatial distribution of the paired closed shell spin up and spin down electrons there is 
no polarisation leading to non-physical effects such as no hyperfine coupling of excess spin- 
density to nuclear spin. On the other hand using an unrestricted approach allows for spin- 
polarization as the spin up and spin down orbitals can have different amplitudes in 
different regions of space. Such unrestricted wavefunctions mix in states of higher S (but 
common Ms). In general the higher this spin-contamination (i.e., how far above S(S+1) the 


expectation value < F >y) the worse the description of the state, and properties, including 
molecular geometries, are less accurate. Thus, one needs to go beyond a single configuration 
to describe these high-spin open-shell states more accurately, and such methods are the only 
starting point for the low-spin open-shell states. For a fuller description see the excellent 
review by Bally et al. which has a chemical discussion of all these aspects.(Bally & Borden, 
2007) Multi-configurational methods, as described in 2.1.3, are therefore essential for a 
balanced treatment of all the states in the oxygen molecule, and this includes when oxygen 
weakly interacts with a sensitizer, the topic covered in the next section. 


4.1 Photosensitisation 


In the discussion of the various aspects of the photochemistry that underpins the pathways 
involved in PDT, we have now arrived at the crucial stage where the reactive singlet oxygen 
species is generated through photosensitisation by the excited state of the irradiated 
photosensitiser molecule. The direct absorption of pure gaseous O2 is very weak due to the 
nature of the electronic transitions, as discussed above, and in order to describe the 
generation of singlet oxygen through sensitisation, the nature of the perturbed O2 molecule 
needs to be considered.(Badger, Wright et al, 1965) As mentioned previously, the 
perturbation is normally viewed as taking place through a collision complex between 
oxygen and the perturbing molecule, in our case the photosensitiser molecule, PS. This 
complex, which we will denote as PS-O,, and a selection of the lowest lying electronic states 
in the complex, is described in figure 11. 


One of the first attempts to describe the mechanisms behind photosensitisation originated 
from observations involving heteroaromatic solvent molecules as perturbers of the O2 
molecule. The fact that oxygen interacts with heteroaromatic organic molecules was first 
inferred in the 1950s when it was shown that samples where oxygen had been dissolved 
gave rise to an absorption band at longer wavelengths than the aromatic compounds alone. 
The first suggestion as to the nature of this interaction, and its precise mechanism, was 
presented in a seminal paper in the 1960s when Mulliken and Tsubomura suggested it to be 
due to absorption by a charge transfer complex (denoted 13(PS*Oz) in Figure 11) originating 
from the interaction between oxygen and the aromatic perturbing molecule, where the 
oxygen has the role as electron acceptor. The interaction in these complexes were considered 
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Fig. 11. Diagram illustrating some states of the PS-O2 complex. As new wavefunctions are 
needed to describe the complex, transitions that were forbidden in the isolated molecules 
can now become allowed. 


to be very small, as the added electron in O, will vacate one of the already occupied mg 
orbitals, which have a very small overlap with the heteroaromatic n-donor 
orbitals.(Tsubomura, 1960) The mechanism of the enhancement of a transition in Oz through 
a collision complex has been suggested to be due to a ‘borrowing’ of intensity from the spin 
allowed transition. The degree of enhancement in these cases will however depend on the 
nature as well as the orientation of the collider, as the degenerate anti-bonding ng and tg y 
will be perturbed differently depending on the nature and the orientation of the valence MOs 
of the colliding sensitiser molecule. In the case of a photosensitiser the size of a porphyrin 
derivative, this is a crucial point in the modelling of the charge transfer event. 


Several more possible mechanisms were suggested through the years, where the main 
theme was the enhanced intersystem crossing of the photosensitiser, where the states are 
mixing either directly or through charge transfer (CT) states of the complex. Kawaoka 
however later suggested that there might not be a need to populate the actual CT states of 
the PS-O2 complex in order to generate singlet oxygen, as its states were mixed to a varying 
degree.(Kawaoka, 1966) It was suggested that CT character will in fact have a negative effect 
on the production of On(a!Ag). This holds true for intramolecular CT within the 
photosensitiser, which facilitates coupling of an excited state to a ground state, as well as for 
intermolecular CT. The latter is normally described as the degree of mixing of the CT states 
with the excited or ground state of the exciplex. (Kawaoka, 1966) If the PS-O2 complex has 
large CT character, radicals are easily formed, which has a negative effect on the singlet O2 
yield. Also, the CT state will facilitate coupling between excited PS-O2 states and the ground 
triplet state 3(!PS-2O2). A further mechanism has also been suggested as it has been shown 
that O2(1A,) can be produced through a direct cooperative absorption by the PS-O2 (32,°) 
complex. (Pagani, Salice et al., 2010; Scurlock & Ogilby, 1988) This highlights the fact that 
many uncertainties still exist with regards to the understanding of the actual mechanism of 
the charge transfer pathway taking in the photosensitisation process. 
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It has become clear that Oz quenches excited states of organic molecules of both 
multiplicities, with varying rates depending on the environment. The quenching of the Sı 
state does not generally contribute to the generation of singlet oxygen, but leads to the 
formation of T; of the sensitiser, through O2 enhanced ISC, or the generation of the sensitiser 
singlet ground state, O2 enhanced IC, with regeneration of the O2(!2Zgt).(Schmidt & 
Schweitzer, 2003) With respect to PDT applications, quenching of the S1 of the sensitiser can 
therefore be looked upon largely as a loss process. 


In terms of modelling the interaction between the oxygen and the sensitizer one is faced 
with several problems. Firstly the weak nature of interaction complex, and the possibility of 
many geometrical conformers contributing in the condensed phase has to be considered. 
Secondly, and the open-shell nature of the collision complex requires the same techniques as 
isolated Oz requires, i.e. multi-configurational methods. Given the often large size of the 
sensitiser this presents many challenges to theory. In fact it is fair to say that accurate multi- 
state, multi-configurational investigations of this system have not been attempted to date. 
Modern methods (for example as discussed in section 2.2.3) can in principle be applied to 
the full interaction complex. These when combined with molecular dynamics simulations, 
and sampling, offer up the possibility that this final crucial step in PDT can be modelled as 
accurately as the other steps and thus hopefully lead to predictive use of modelling in the 
design of new PDT drugs. 


5. Conclusions 


We have introduced Photodynamic Therapy in terms of the photochemical pathways taking 
place within a target biological system. We presented an outline of the challenges associated 
with the modelling photosensitiser molecules, focusing on the issues related to electron 
correlation problem, in terms of methods using the electronic wavefunction as well as 
methods based on the electronic density of the system. Further we have highlighted the 
computational advances that have made it possible to accurately describe both linear and 
non-linear absorption processes in molecules the size of common photosensitisers. An 
accurate description of these processes can be crucial in the development of photosensitisers 
with absorption maximum within the optical window of tissue penetration. The modelling 
of ultrafast processes taking place from an excited state of a molecule is challenging and we 
have summarised the issues associated with the description of these processes in large 
molecules. Finally, we have discussed the difficulty in an accurate theoretical description of 
the energy transfer step involving molecular oxygen present in the target cell, mainly in 
terms of the challenges associated with the modelling of open-shell systems. 


Despite several drawbacks, including poor light penetration of cancer tissue and poor 
targeting of the photosensitiser systems, PDT remains a very attractive goal, and research 
continues through many avenues. The search for the ideal photosensitiser, which will make 
PDT a realised anti-cancer treatment for general use is very much under way and the input 
that can be gained from computational studies cannot be ignored, as it has the potential to 
increase the understanding of the entire photochemical pathways involved. Computational 
techniques have now matured so that molecules the size of standard photosensitisers can be 
investigated routinely for their absorption properties. This is also true for non-linear 
absorption properties, which means that this approach to PDT can be developed and hence 
potentially improve the scope of the technique, and be an important tool both when 
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rationalising the existing classes of photosensitizers and the design of future 
photosensitisers with desired molecular properties. 


Even though significant challenges clearly still remain, the insight a theoretical input can 
give to the understanding of the fundamental photochemical mechanisms has potential to 
greatly accelerate the development of PDT as an everyday technique. 
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